### Fig. 5. Composition of NCS constraints in ten countries with highest NCS biophysical mitigation304potential
library(tidyverse)
library(dplyr)
library(ggplot2)

nbase<-read.csv("NBaseData_v1.csv")
NCS_potential <- read.csv('NCS_potential Fig 5.csv')

constraint_category <- nbase %>%
  distinct(Constraint_short, Constraint_category) %>%
  group_by(Constraint_category) %>%
  summarise(total_constraint_count = n()) %>%
  ungroup()


NCS_potential$width <- NCS_potential$NCS_potential/600

# Summarize constraints by country and category
country_constraint_category <- nbase %>%
  distinct(Country_code, Constraint_short, Constraint_category) %>%
  group_by(Country_code, Constraint_category) %>%
  summarise(category_count = n()) %>%
  ungroup()

# Merge with NCS potential data
country_constraint_category <- merge(country_constraint_category, NCS_potential, by = "Country_code", all.y = TRUE)

# Calculate proportions within each country
country_constraint_category <- country_constraint_category %>%
  group_by(Country_code) %>%
  mutate(proportion = category_count / sum(category_count))

# Merge the total category count for each constraint category to show in the legend
country_constraint_category <- country_constraint_category %>%
  left_join(constraint_category, by = "Constraint_category")

# Create a label with both the category name and total number of constraints
country_constraint_category <- country_constraint_category %>%
  mutate(label = paste0(Constraint_category, " (n=", total_constraint_count, ")"))

country_constraint_category$potential_width <- country_constraint_category$NCS_potential/600

country_constraint_category_data <- country_constraint_category %>%
  filter(!is.na(Constraint_category))

country_constraint_category_nodata <- country_constraint_category %>%
  filter(is.na(Constraint_category)) %>%
  mutate(proportion = 1,
         Constraint_category = "No Constraint data",
         label = "No Constraint data")

country_constraint_category_rbind <- bind_rows(country_constraint_category_data, country_constraint_category_nodata)

# Plot the horizontal stacked bar chart with varying widths (mitigation potential)
NCS_potential_label <- country_constraint_category_rbind %>%
  distinct(NCS_potential) %>%
  arrange(NCS_potential)

country_constraint_category_rbind$Source_label <- paste0("(N = ", as.character(country_constraint_category_rbind$Source), ")")

plot <- ggplot(country_constraint_category_rbind, aes(y = reorder(Country_code, NCS_potential), x = proportion, fill = label)) +   
  
  # Bar for countries with constraint data
  geom_bar(stat = "identity", width = country_constraint_category_rbind$potential_width, aes(fill = label), color = "black") +
  
  # Bar for countries without constraint data (No data), hide "No data" from the legend
  geom_bar(data = country_constraint_category_nodata, aes(y = reorder(Country_code, NCS_potential), x = proportion), 
           stat = "identity", width = country_constraint_category_nodata$potential_width, fill = "#dddddd", color = "black", show.legend = FALSE) +
  
  # Add the count of constraint category to the middle of each stack
  geom_text(data = country_constraint_category_data, aes(label = category_count), 
            position = position_stack(vjust = 0.5), size = 5 / .pt, color = "black", family = "Helvetica") +  
  
  # Add the "No data" text inside the white stacks for countries without constraint data
  geom_text(data = country_constraint_category_nodata, aes(x = 0.5, label = "No Constraint data"), 
            size = 5 / .pt, color = "black", vjust = 0.5, family = "Helvetica") +
  
  # Add the country name to the right of each stack
  geom_text(aes(x = 1.05, label = Country_code), hjust = 0, size = 5 / .pt, nudge_x = -0.04, nudge_y = 0.4, family = "Helvetica") +  # Move country name slightly to the right of the bar
  
  geom_text(aes(x = 1.05, label = Source_label), hjust = 0, size = 5 / .pt, nudge_x = -0.04, family = "Helvetica") +
  
  # Customize the Y-axis to display NCS potential instead of country names
  scale_y_discrete(labels = NCS_potential_label$NCS_potential) +
  
  # Add NCS Potential (MtCO2) to the Y-axis label with proper subscript
  labs(x = "Constraints Composition", y = expression("NCS Mitigation Potential (MtCO"[2]*"/yr)"), fill = "Constraint\n Category") + 
  
  # Manual color scale for the constraint categories
  scale_fill_manual(values = c("Finance (n=6)" = "#77aadd", 
                               "Knowledge (n=7)" = "#99ddff",
                               "Material Inputs (n=4)" = "#44bb99",
                               "Social and Behavioral (n=9)" = "#bbcc33",
                               "Government and Organizations (n=5)" = "#aaaa00",
                               "Markets (n=7)" = "#eedd88",
                               "Negative side effects of NCS (n=1)" = "#ee8866",
                               "Rules and Laws (n=7)" = "#ffaabb")) +
  
  # Theme settings to remove grey background, remove grid lines, and add axis lines
  theme_minimal() +
  theme(
    text = element_text(family = "Helvetica", size = 5),
    panel.background = element_blank(),  # Remove the grey background
    panel.grid = element_blank(),        # Remove all grid lines
    axis.line = element_line(linewidth = 0.1 * .pt, colour = "black"),          # Add visible X and Y axis lines
    axis.text.x = element_blank(),       # Remove X-axis labels
    axis.ticks.x = element_blank(),      # Remove X-axis tick marks
    axis.ticks.y = element_blank(),      # Remove Y-axis tick marks but keep the labels
    legend.position = "top",
    legend.text = element_text(size = 5), # Change font size of text in legend only
    legend.title = element_text(size = 5),
    axis.title.x = element_text(size = 5), # Increase font size of X-axis label
    axis.title.y = element_text(size = 5)
  )

ggsave("Fig 5.pdf", plot = plot, width = 180, height = 135, unit = "mm", dpi = 300)

